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Abstract 

This paper describes algorithms to deal with nested symbolic sums over combinations of har- 
monic series, binomial coefficients and denominators. In addition it treats Mellin transforms 
and the inverse Mellin transformation for functions that are encountered in Feynman diagram 
calculations. Together with results for the values of the higher harmonic series at infinity the 
presented algorithms can be used for the symbolic evaluation of whole classes of integrals that 
were thus far intractable. Also many of the sums that had to be evaluated seem to involve new 
results. Most of the algorithms have been programmed in the language of FORM. The resulting 
set of procedures is called SUMMER. 
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1 Introduction 



The computation of Feynman diagrams has confronted physicists with classes of integrals that are 
usually hard to be evaluated, both analytically and numerically. Also the newer techniques applied 
in the more popular computer algebra packages do not offer much relief. Therefore it is good to 
occasionally study some alternative methods to come to a result. In the case of the computation 
of structure functions in deep inelastic scattering one is often interested in their Mellin moments. 
Each individual moment can be computed directly in ways that are much easier than computing 
the whole structure function and taking its moments afterwards. There exist however also instances 
in the literature in which all moments were evaluated in a symbolic way ||] ||] || [||. Once all 
positive even moments are known, one can reconstruct the complete structure functions. Hence such 
calculations contain the full information and are in principle as valuable as the direct evaluation 
of the complete integrals. In these calculations the integrals become much simpler at the cost of 
having to do a number of symbolic sums over harmonic series. The draw-back of the method is 
that although much effort has been put in improving techniques of integration over the past years, 
very little is known about these classes of sums. A short introduction is given for instance in 
ref ||. In addition such calculations are of a nature that one needs to do them usually by means 
of a computer algebra program. This means that when algorithms are developed, they should be 
suitable for implementation in the language of such programs. 

This paper describes a framework in which such calculations can be done. As such it gives a 
consistent notation that is suited for a computer program. It shows a number of sums that can 
be handled to any level of complexity and describes an implementation of them in the language 
of the program FORM Q. Then the formalism is applied to the problem of Mellin transforms of 
a class of functions that traditionally occurs in the calculation of Feynman diagrams. This in its 
turn needs harmonic series in infinity and hence there is a section on this special case. Next the 
problem of the inverse Mellin transform is dealt with. With the results of the series at infinity one 
can suddenly evaluate a whole class of integrals symbolicaly. This is explained in the next section 
where some examples are given. 

The paper is finished with a number of appendices. They describe the details of some of the 
algorithms and their implementation. Additionally there is an appendix with lists of symbolic sums 
that are not directly treated by the 'general' algorithms. These sums were obtained during various 
phases of the project and many of them do not seem to occur in the literature. 

2 Notations 

The notation that is used for the various functions and series in this paper is closely related to how 
useful it can be for a computer program. This notation stays as closely as possible to existing ones. 
The harmonic series is defined by 




(1) 




(2) 



in which m > 0. One can define higher harmonic series by 




(3) 
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S-m,n,-, jp {n) = Y. { -^S h ,., jp (i) (4) 
i=i 1 

with the same conditions on m. The m and the ji are referred to as the indices of the harmonic 
series. Hence 

i=l j=l •' fc=l 

In the literature the alternating sums are usually indicated by a bar over the index. The advantage 
of this notation is that it can be extended easily for use in a computer algebra program, eg.: 

Sh,-,i m {n) -» S(R(il, . . . ,im),n). 

Such objects can be easily manipulated in the more modern versions of the program FORM. 

The argument of a harmonic series which has only positive indices can be doubled with the 
formula: 

S mi ,-,mM = E2 mi+ - +mp - p ^ 1 ,.., ±Jp (2n) (6) 
± 

in which the sum is over all 2 P combinations of + and — signs. 

The weight of a harmonic series is defined as the sum of the absolute values of its indices. 

m 

W(S jl ,., jm (n)) = ^\M (7) 
i=i 

For any positive weight w there are 2 x 'i w ~ 1 linearly independent harmonic series. The fact 
that for each next weight there are three times as many can be seen easily: One can extend the 
series of the previous weight either by putting an extra index 1 or —1 in front, or by raising the 
absolute value of the first index by one. 

The set of all harmonic series with the same weight is called the 'natural' basis for that weight. 

The extended weight of the compound object of a series and denominators is the weight of the 
series plus the number of powers of denominators that are identical to the argument of the series. 
Hence S\^^{n)/n A has the extended weight 13. 

The total weight of a term is the sum of all extended weights of all the series in that term. 
Hence the total weight of 82,3(11) S-2(m) is 7. 

The value for an index is reserved for an application that is typical for computers. If the 
results of a given weight need to be tabulated, the above notation would require a table in which 
the number of indices is not fixed. This can be remedied by a modified notation which is only 
used in specific stages of the program. An index that is zero which is followed by an index that is 
nonzero indicates that one should be added to the absolute value of the nonzero index. Hence: 



50,1,0,0,-1,1 ( n ) = s 2,-3,i(n) (8) 

This way one can express all series of weight w into functions with w indices of which the first w—1 
can take the values 1, and —1, while the last one can take the values 1 and —1. 

Consider the following identity which can be obtained by exchanging the order of summation: 

Sj, k (n) + S kJ (n) = Sj(n)S k (n) + S jhk {n) (9) 

in which the pseudo addition operator & adds the absolute values and gives the result a positive 
value if j and k have the same sign and otherwise the result will have a negative value. One could 



2 



select a basis in which one keeps products of harmonic series with as simple a weight as possible. 
The above equation would indicate that in that case one of the left terms should be excluded from 
the basis in favor of the first term on the right hand side. Although the choice of which of the 
higher harmonic series to keep and which to drop in favor of the product terms is not unique, there 
are cases in which such a basis is to be preferred. In particular when n — > oo one can choose a 
basis in which all divergent objects are expressed as powers of 5*1(00) multiplied by finite harmonic 
series. In general however the summation formulae are much simpler in the natural basis in which 
each element is a single higher harmonic series. This can be seen rather easily when looking at the 
sum: 

■A Sj(i)S k (i) _ A S jt k(i) + Sk,j(i) - S jkk (i) 

1=1 j=l 

= S m j tk (n) + S mj k,j(ri) - S m jk k (ri) (10) 

In order to be able to do the sum one has to convert to the natural basis anyway. After the summa- 
tion one would have to convert back to whatever other basis one happens to work with. Appendix 
A presents an algorithm by which combinations of harmonic series with the same argument can be 
expressed in the natural basis. 

Additionally there can be denominators containing the summation parameter. There arises 
immediately a problem when there is more than one denominator. Traditionally one can split the 
fractions with 

11 1 .1 1 . 



i + a i + b b — a i + a i + b 
Unfortunately this formula is not correct when a = b. Because often there will be nested sums and 
sums with symbolic parameters, a and b can be functions of summation or other parameters and 
hence it will not be obvious when a = b. In FORM this can be repaired in principle with one of 
the special functions: 

11 1 ,,.,1,1 1 , 

i+ai+b (i+a) z b — a 1 + a i + b 

Here 5 a ^ = 1 when a = b and zero otherwise. In the language of FORM it is represented by the 
built in object delta_. Unfortunately this form of the partial fractioning is not very useful, because 
it still evaluates into terms involving l/(b-a) in which a can be equal to b. Hence an even more 
complicated form is needed: 

! ' \ b7 r^ + (9(a-b-l)+e(b-a-l))-^—(-^—-^-) (13) 



i+ai+b ' (i+a) 2 b — a i + a i + b 

in which one has to assume that a and b only take integer values. The function 6(x) (in FORM 
theta\_(x)) is zero when x is negative and one when x is zero or positive. These ^-functions fulfill 
the role of conditions like a > b+1 plus b > a+1 and are worked out first. Hence this should not 
be read as 0/0 for the case that a = b. The complete and proper equation would involve a new 
function: 

— ^7 = 6 ab ^-^ + (9'(a-b) + 6'(b-a))-!—(— (14) 

i+ai+b ' (i+a) 2 y y > v ^5_ a ^ + a i + > v > 

with 6'(x) (in FORM thetap\_ (x) ) is one when x > and zero when x < 0. Actually 8'(x) = 
l — 6(—x), but this cannot be used for the same reason that the equation (|l~2"l) could not be used. 
Because it is rather complicated to manipulate both the functions 9(x) and 9'(x) simultaneously, and 
because one has almost always integer values of the parameter, the computer program uses mostly 
the formula (13) which assumes the integer values. It should be clear that much attention should 
be given to theta and delta functions, their combinations and their interactions with summations. 
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3 Synchronization 



When one has to do sums over a combination of objects one of the problems is that such objects 
do not always have identical arguments. If this is the case one would have to program many more 
sums than often is necessary. Whenever it is possible one should 'synchronize' the arguments. 
This means that one tries to make the arguments of the various harmonic series, the denominators 
and the factorials equal to each other. This can be illustrated with one harmonic series and one 
denominator: 

t=i i=i i=i 

s_m i _i_ (15) 

1=1 1=1 1=1 

In this equation and the sequel it is assumed that the left most index is positive. If it is negative 
there will be the extra (— l) 1 and one has to be more careful with the signs of the terms, but the 
principle is always the same. 

Of course, when the difference between some arguments is symbolic like in S\(i+k)/i, such tricks 
do not work, but for differences that are integer constants one can define a scheme that converges. 
Let m be a positive integer constant in the remaining part of this section. In that case one can 
write: 

Sj,r 1 ,-,r a ( i + m ) = ff?,ri,-,r„(* + ra-l) | S ri ,...,r a (.i + m) 

i i i(i+m,y 

The partial fractioning of the denominators in the last term results in terms that have only a power 
of l/[i+m) and one term which has a factor This last term however has a simpler harmonic 
series in the numerator. Hence this relation defines a recursion that terminates. Similarly one can 
write: 

Sj,n,-,r a (i) = S jtru ... trs (i + 1) _ S ri ,...,r a (.i + l) ^ 

i + m i + m (i + m)(i + \)i 

and partial fractioning results again in terms in which the arguments either are the same, or closer 
to each other, or the harmonic series has become simpler. 
Next is the interaction between two harmonic series: 

Sj,r 1 ,---,r s {'i)S Ply ...,p q {i+m) = Sj irij ... jrs (i+l)Sp lj ... jPq (i+m) 

-S ri ,-,r s {i+l)S pu ..., Pq (i+m) * (18) 

This relation defines, in combination with the previous two equations, also a proper recursion. In 
the last term one can synchronize the argument of the second harmonic series with that of the 
denominator, giving (potentially many) terms with either l/(i+m) or an argument that is closer 
to i+1. In all cases the arguments are at least one closer to each other. In addition some of the 
harmonic series have become simpler. 

Once two harmonic series have the same arguments this product can be rewritten into the basis 
of single higher harmonic series (see appendix A) . Hence products of more than two harmonic series 
with different arguments can be dealt with successively. 
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At this point there can still be factorials. The beginning is easy: 



i (i+m)\ i (i+m) (i+m— 1)! 

1111 

= (19) 

mi (i+m— 1)! m(i+m)\ 



and 

1 i+l 



(i+m) i! (i+m) (i+l)! 

= _J ^Zl (20) 

(i+l)! (i+m) (i + l)! 1 ' 

Because of these two equations one can also synchronize combinations of harmonic series and 
factorials. 

The problem is that usually one cannot do very much with the product of two factorials. 
This means that if one has more than one factorial, one may be left with factorials with different 
arguments. 

Another problem exists with arguments of the type i versus arguments of the type n—i. These 
can of course not be synchronized completely, but if n is the upper limit of the summation over i, 
one can try to make a synchronization that excludes other nonsymbolic constants. This is slightly 
more complicated than what was done before: 

— S jtri ,., rs (i+1) = —Sj,r u -,r s (i) + ( n _ i )( i + 1 )j S '- 1 ,.,r.(* + l) (21) 

Partial fractioning of the last term will leave something simpler. Similarly there is: 

—S jjri: ... :rs (i-1) = —.S j , ri ,.. tr ,{i) - (n _. )( ._ 1)j 5 ri) .., J . s (i-l) (22) 

For two 5-functions one can write: 

Sj,r 1 ,--,r 3 ( n —i)Sk, Pl ,---, Ps (i+^) = Sj^ ly ... yrs (n — i)Si t .p 1 ^....p s (i) 

+ (j_)_l)A:^' ri '-' rs ( n_ ^'^P 1 '-'P s ^ + 1 ) 
= ^j,ri,---,r a ( n ~^)^k,pi,---,p s (i) 

+ (j+l)fc^'' ri ''''' r ^ n_i_1 ^ pl '''''P s ^ + 1 -' 

+ ( i+ l)A ; 1 ( n _ i )i ^i^^( n - i ) 5, W^^( i + 1 ) ( 23 ) 

Again partial fractioning of the last term leads to a simpler object. One can derive equivalent 
relations for combinations involving factorials. In this case also pairs of factorials can be dealt 
with: 

1 1 1,11, 

+ 



(n—i)! (i+l)! (n — i — 1)! i! n+1 n—i i+l 

^+T^(n-i)! i! + (n-i-1)! (i + 1)!^ ^ 

All the above relations can be combined into one recursion that leaves all S-functions, all denomina- 
tors and at least one factorial properly synchronized. Additionally one has a proper adjustment to 
the boundaries of the summation, and therefore the factorials can often be combined into binomial 
coefficients. 
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4 Mellin Transforms 



The Mellin operator M is defined by 



and the operator M + by 



with 



M{f(x)) = f 1 dx x m f(x) (25) 
Jo 

M + (f(x))= f 1 dxx m -^- (26) 
Jo (1-^)+ 

when /(l) is finite. When there is a power of ln(l— x) present it becomes 

M+((ln(l-x)) fc /(x)) = (\xx m {H )- x)) * f{x) (28) 

Jo \l- x )+ 

1 dx(f{z) - f(l)) ^- x))k (29) 
l—x 

These are the traditional operations. In the literature one often defines the transform shifted over 
one as in 

M(/(x)) = f dx x N - l f(x) (30) 
J 



In the context of this paper the notation will be the one of equation (25). The 'Mellin parameter' 
is given in that case in lower case variables. Hence the translation to the shifted notation should 
be of the nature n — > N— 1. 

For Mellin transforms of formulas resulting from Feynman diagrams one has to consider the 
transforms of functions that are combinations of 1/(1— %)+, l/(l+a?), hi(x), ln(l+x), ln(l— x), powers 
of these logarithms, and various polylogarithms of which the arguments are rational functions of 
x. Powers of x just change the moment of the function. Hence they do not have to be considered. 
Additionally one can always assume that either 1/(1 — 2;)+ or 1/(1 + x) is present, because the 
functions without such a term can be written as two functions in the class that is being considered: 

1 x , 

1 = + (31) 

l+x 1+x v ' 

The algorithm that obtains the Mellin transform of any combination of such functions is rather 
direct. Consider the following steps: 

1. If there is a power of 1/(1 — x) or l/(l+x), replace it by a sum according to the formulas 

= s> f < 32 > 

CO 

= (-ir£(-i)v (33) 



1 + X 

I Y/fc 

2. If the function to be transformed contains powers of ln(l— x), split it into its powers of In (l—x) 
and F(x) which represents the rest and has a finite value at x = 1. Then one writes 



f 1 dx x m ln p (l-x) F(x) = C dx x m ln p (l-x) (F(x) - F(l 
Jo Jo 



Jo 



)) 



dx x m ln p (l-x) (34) 
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3. The Mellin transform of just a power of ln(l — x) can be replaced immediately using the 
formula 

/ dxx m \n p {l-x) = K ' P 5i-i(m+l) (35) 
jo m+1 ' ' 

in which the S- function has p indices that are all 1. This avoids divergence problems during 
the next step. Similarly one can apply: 

/ dx x m \vP{x) = . K ,P (36) 
Jo V ' (m+1)^ 1 v ' 

when there is only a power of ln(x) left, but this step is not essential; it only makes the 
algorithm a bit faster. Due to the powers of x there will be no divergence problems near 
x = 0. 

4. Do a partial integration on the powers of x. Because of the second step, the values at x = 
and x = 1 never present any problems. 

5. If there is only a power of x left one can integrate and the integration phase is finished. 
Otherwise one should repeat the previous steps until all functions have been broken down. 
Note that for this to work all functions have to break down properly. Hence one cannot use 
fractional powers of the functions involved. 

6. At this point the terms may contain nested sums, either to a finite upper limit or to infinity. 
These sums do not present any complications once products of two ^-functions with identical 
arguments can be combined into elements of the natural basis (see appendix A). 

The main complication in the above algorithm is the treatment of the infinities that may arise in 
the summations. Many of the terms develop a divergence. These are all of a rather soft nature 
and hence their regularization is relatively easy. All divergences in the sums are of a logarithmic 
nature and hence, if one considers the sum to go to a rather large integer L, the divergent sums 
behave like powers of InL up to terms of order 1/L[]. Because all transforms should be finite the 
terms in InL should cancel. After that one can safely take the limit L — > oo. Taking this all in 
consideration, all sums that contain a divergence can be rewritten into powers of one single basic 
divergent sum (Si(oo)) and finite terms. After that there are no more problems of this nature. 

The result of the above algorithm is an expression with many harmonic series of which the 
argument is a function of m and others of which the argument is infinity. These last sums are 
treated in the next section. 



5 Values at infinity 

In the previous section the results of the Mellin transforms were harmonic series in the Mellin 
parameter m and harmonic series at infinity. In order to solve the problem completely one has to 
find the values for these series at infinity. After all they represent finite numbers and the number of 
series is much larger than the number of transcendental numbers that occur once they are evaluated. 



The sums to be considered are related to the Euler-Zagier sums |9[ |10[| which are defined as 



" a- 

C(s 1 ,...,s k ;a 1 ,...,a k )= ^ 11^7" ( 37 ) 



Tl ■ 

rij>nj + i>0 j=l J 



1 In principle there is also an Euler constant, but when the logarithms cancel, also the Euler constants cancel and 
hence they are not considered here 
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These sums are however not identical to the 5-functions at infinity because for them the sum is 



* »M- E <*) 

nj>n j+1 >l j=l Tlj 



The notation [ ] sl <o indicates that this part is present only when s\ < 0. Here a method is 



presented to evaluate these sums that is completely different from the one in reference [11|. 

The first step in the evaluation of the sums is to express the sums as much as possible in terms 
of products of harmonic series with a lower weight. This can be done up to a point. One will always 
need a number of series with the weight one is considering. This step is basically the inverse of the 
algorithm of appendix A. It is harder to be implemented in a deterministic way, because the choice 
of the basis is not unique. But this can be solved in a different way as will be seen below. 

Next there are two types of extra identities one can consider. The first set comes from looking 
at the series with only positive indices and applying the doubling formula (||) to it. For all the 
series that are finite it makes no difference whether the argument is infinity or two times infinity. If 
the selected basis is such that all divergences are powers of <Si(oo) one only has to make the extra 
adjustment Si(2oo) — > 61(00) + ln(2). This gives a number of extra equations that correspond to 
new relations between the series. Unfortunately this does not give enough relations, but some are 
interesting in their own right. For instance 

5 m (oo) = 2 m - 1 (5 m (2oo) + 5_ m (2oo)) (39) 

gives immediately the well known relations 

S_ m (oo) = -(1 - 2 1 - m )5 m (oo) m > 1 (40) 

5_i(oo) = -ln(2) (41) 

The more powerful consideration however is the following: Suppose one is summing over a 
square grid of size n x n. Under what conditions is the sum over the upper right diagonal half 
of the square {i\ + 12 > n) zero in the limit n — > 00? If this sum is zero, the product over two 
individual sums can be replaced by a sum over the lower left diagonal triangle (i\ + 22 < n). This 
leads to the following theorem: 

Theorem: When not both mi = 1 and k\ = 1 the following identity holds: 



S mi ... tmv (oc)S kl „. <kq (oc) = \^J2 S mir-,m P {n-i)S k2j ... tkq (i) 1 J* l L<0 (42) 

i=i 



The proof is rather trivial, considering that all m^ and k{ are integers and that alternate series with 
(— 1)* actually converge one power of i better than they seem to at first sight. This can be seen 
when the terms are grouped in pairs. The sums can be estimated by integrals and the numerators 
can only give powers of logarithms. Hence the presence of at least three powers of denominators 
(excluding m\ = k\ = 1) will make the limit go to zero. 

The sum can be readily worked out with the algorithm described in appendix C. 

Assume that all sums up to weight n have been determined. The complete algorithm for weight 
n + 1 is now: 

1. Construct all pairs of S'-functions for which the sum of the two weights is n + 1. 

2. Each pair is used to construct two equations (unless both S'-functions have their first index 
equal to one in which case the second equation that would have been based on the above 
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theorem is not made). The first equation is made by taking S^S^ — S^S^ and applying 
the routine (see appendix A) that converts the 5-functions to the basis to the first pair. 
These are the 'shuffle algebra' relations. The second set of equations is created by taking 
gC 1 ) <>( 2 ) _ gWgl 2 ) and then applying the formula of the theorem to the first pair. After this 
the routine of appendix C is applied. 

3. Substitute the values for the lower S- functions. 

4. Eliminate now the 'unknown' S'-functions of weight n + 1 as much as possible as if one is 
solving a linear set of equations (which is what it is). Apply the same set of substitutions 
that will eliminate the equations to the series that need to be evaluated. 

5. Inspect the result and see which sums should be considered as new independent variables 
because they were not eliminated. If one insists on a given sum to be among the variable(s) 
not to be eliminated one can substitute it by a different variable before the elimination 
procedure. 

It is not so difficult to construct a program in the language of FORM that can execute this procedure 
all the way to S'-functions of weight 7. Such a program takes just a few hours (< 6 without special 
optimizations) on a Pentium-II-300 processor. When a series diverges one uses the basic divergence 
5*1(00) as if it were a regular variable. This presents no problems. 

The variables that one needs at the different weights are: S*i(oo), ln(2), £2, C3> LUQ), (5, 
1,15(2), ^6(q)j *5'-5,-i( 00 )) C7j ^7(2), 5'-5,i,i(oo), 5 , 5 i _i j _i(oo). The choice of the S'-functions that 
remain is not unique. Here the selection is such that they contain as few indices as possible and 
are as convergent as possible. Numerical values for these quantities can be obtained by standard 
techniques. 

LU(l/2) = 0.51747906167389938633 
Li 5 (l/2) = 0.50840057924226870746 
Li 6 (l/2) = 0.50409539780398855069 
Li 7 (l/2) = 0.50201456332470849457 
S-5-i(oo) = 0.98744142640329971377 
-£-6,1,1(00) = 0.95296007575629860341 

5*5 _i_i (00) = 1.02912126296432453422 (43) 



It should be noted however that according to the work by Broadhurst and Kreimer [12] most of these 
constants should not appear in the computation of massless Feynman diagrams. The first non-zeta 
constant should be 65,3(00) which is an object of weight 8. This indicates that in x-space the 
functions can only occur in such combinations that these constants cancel in Mellin space. Hence 
one may not need to know their values for many applications. In the case of massive Feynman 
diagrams the situation is different. The constant Li4(l/2) does occur in the three loop corrections 
to the g—2 of the electron [13]. 



The results of the runs up to weight 7 have been tabulated and put in the FORM program. The 
main problem in making the tables is that the objects with identical weights may have different 
numbers of indices. Hence the notation of indices that are either —1, 1 or of equation (||) is used 
for the tables. The conversion to and from this notation is rather simple. 
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6 Inverse Mellin Transforms 



If one can obtain a result in Mellin space (as a function of n) in principle it is possible to convert 
to the function in x-space. This is however a rather complicated operation. There exists some 
literature about it [|j Q but it remains rather difficult. Also considering it as some type of Laplace 
transform does not give much relief || . In many cases one can employ a different strategy. Given a 
result in Mellin space with a set of series, one can try to find a set of functions in x-space for which 
the Mellin transforms span the space of the functions in Mellin space. After that one only has to 
solve a set of linear equations to make the inverse transform. In the case of two loop moments of 
structure functions in deep inelastic scattering, the results in Mellin space are just 5-functions of 
weight 4. Because the whole space of such S'-functions is 54 dimensional (a basis has 54 elements) 
one has to find 54 functions in x-space that map into the Mellin space in a linearly independent way. 
This does not present too many problems. One should of course note that this method depends on 
having routines to do the Mellin transforms automatically. 

For higher weights it may not be so easy to find a complete set of functions in x-space. This can 
be illustrated by a simple calculation. To obtain a complete set of functions in x-space for which 
the Mellin transforms cover the natural basis of weight w one needs 2 x S""^ 1 functions in x-space. 
Because this number can be divided by two (the relevant functions are of the types /(x)/(l ± x)) 
only S"^ 1 functions have to be considered. A number of these can be constructed by taking products 
of functions that contribute to lower weights. That leaves a number of functions that are new at 
the given weight. This number increases rapidly with the weight. They are 3,8,18,48,116 for 
the weights 3,4,5,6,7 respectively. Hence one has to come up with a rather large number of new 
functions when the weight becomes large. Fortunately there is a method that will work provided 
only a numerical answer is needed for any value of x. 

Assume that for a given weight w all necessary functions in x-space are known. Assume also 
that the Mellin transform of some F is given by 

£** m = igffi (44) 

in which rft represents any allowable series of the type mi, • • • , m q and p > 0. For this function F 
one has 

/^^S = (-l) n ^,Tn>W-(-l) n ^(00) (45) 

Jo dX (l-x) = S P^(°°^ ~ S P,nt( n } ~ S i(oo)F(l) (46) 

In the second expression one can see that -F(l) will be nonzero when p = 1 and zero otherwise. 
This is needed to keep the expression finite. It is assumed here that F(x) does not contain a factor 
ln(l — x) or that if it does the other components of F still make that F(l) = 0. If this is not the 
case there will be more complicated sums of the type of appendix [D] and the right hand side will 
have more terms to cancel the divergences that are due to S p ^(oo) having more than one power 
of 51(00). Rather than using the sums of appendix |D| one can also use the algorithms of section ^| 
to break down the function F completely. 

Considering that a knowledge of all odd or all even moments is sufficient to reconstruct F the 
presence of (— l) n should not be a problem in the end. It does not lead to a doubling of the necessary 
functions -even moments in terms of N correspond to odd moments in terms of n— . One should 
also observe now that the functions F(x)/(l+x) and i ? (x)/(l— x) are related to the inverse Mellin 
transforms of 5^ ^(n) and S p ^(n) respectively. Assume now that the S p -^{n) are of weight w. 
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How does one construct the inverse Mellin transforms of functions of weight w + 1? For this one 
should have a look at the functions 

F+(x) = [ X dx^± (47) 
JO 1+x 

rx F(r) 

F~(x) = / dx^-t- (48) 
Jo 1— x 

F°(x) = ( l dx^- (49) 

Jx X 

For these functions one can derive readily by means of partial integration 



dx x n F~(x) = —^Sp^in + l) (51) 



r i _1 ytrx 1 

/ dxx n F+(x) = - { -J—S^(n+l) + — ((-1)^-1)5^(00) (50) 
Jo n+l v ' n+l F ' 

1 

o n+l 

rl 1 

/ dx x n F°(x) = —, — -5=*(n+l) (52) 

Jo (n+l)P +1 my ' y > 

With the aid of equations (|Iq) and ( f46[) one derives now the relations 

^""^ = -(-l) n Si rf> ?^) + (-1)^ 

+(-l) n S_ P)? ^(oo)(Si(n) - SL x (n) + £4(00)) (53) 

^""^ = (-l) n (5_ liPi ^(n) - &- 1)P ,^(oo)) (54) 
1 F°fa;) 

dxxn T^ = (- 1 ) n (-^Pfi),^( n ) + 5 -(^-i),^( 00 )) ( 55 ) 

f 1 x"F+(x) 

Jo (1— x)+ ~~ "-l.-jj^W ~~ d-it-p^K 00 ) 

+ s ^p,wt( 00 )(S-i( 00 ) ~ S-i(n) + Si(n)) (56) 

i x n F~(x) 

dx m_ x \ = - S i, P M( n ) + 5 i, P ,^(°°) - Si(oo)S p tf(po) (57) 
i x n F°(x) 

^ ( 1 _ a ) = s (pn),T7t( n ) ~ S (pn),^(°°) ( 58 ) 

In these expressions is assumed that F(x) contains no factors ln(l — x). In that case it is not not 
difficult to see that all divergences cancel. When there are factors ln(l— x) the expressions become 
a bit more complicated in the constant terms in order to obtain a complete cancellation of the 
divergences. The first terms of the right hand side expressions form indeed a complete set of S- 
functions of weight w+1 when all possible values of p and all possible S-functions in equation 
are considered. Because all other terms in the right hand side expressions are of a lower weight in 
terms of the argument n, their inverse Mellin transforms are supposed to be known and hence all 
inverse Mellin transforms of weight w+1 can be constructed. If the integrals in the definitions of F + , 
F~ and F° cannot be solved analytically, one can still obtain their values numerically by standard 
integration techniques. If one has to go more than one weight beyond what is analytically possible, 
one obtains multiple integrals. Many of these can of course be simplified by partial integrations as 
can be seen in the following formula: 

F++(x) = f [ X ^-F{x) (59) 
V 1 Jo 1+xJo 1+x v 1 y ' 
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= in ( i + x) r *l f{x) - r n*) dx m 

Jo 1+x Jo 1+x 

At this point it seems best to give some examples. First look at the constant function in Mellin 
space. It is the only function with weight zero and its inverse Mellin transform is 5(1 — x). Here 
5(x) is the Dirac delta function. Hence the inverse Mellin transforms for functions with weight one 
are: 

(-lfS-i(n) - -1- + (_!)» in(2)*(l-x) (61) 

Si(n) - -J- (62) 
1 — x 

The factor (— l) n in the right hand side indicates that the reconstruction from the even moments 
different from the reconstruction from the odd moments. This means that if the moments are 
obtained for even values of N (which means odd values for n) one should treat the terms in 5^>(n) 
differently from the terms in Sj^(n + 1). 

Next are the functions with weight 2. The only function with weight one that can occur in 
equation (|44|) is l/(n+l) and its inverse Mellin transform is given by F(x) = 1. From this one can 
construct F + (x) = ln(l+x), F~(x) = — ln(l— x) and F°(x) = ln(x). One can now work out the 
equations ( |53| - |58| ) to obtain the inverse Mellin transforms for the weight two functions. 

For the weight three functions one obtains dilogarithms with the arguments x, —x and (l+x)/2 
as new objects. For the weight four functions the functions F (x) and F°(x) can have trilogs with 
the arguments x, —x, (l+x)/2, l/(l+x), 1—x, 2x/(l+x), (1— x)/(l+x) and— (1— x)/(l+x). Of course 
one may choose a different representation in which the function S\^(x) plays a role (see references 



14J and 

There is one more important observation to be made. The expressions (53-|58) have just a single 



S"-function of weight w+1 in the right hand side. This means that one can obtain the inverse Mellin 
transforms of the various S- functions without having to solve sets of equations. One only has to 
move terms from the right hand side and put their inverse Mellin transform (which is much simpler) 
into the various F-functions. This can be done systematically and it can be checked by the Mellin 
transformation program. The approach of looking for which functions can occur and then making 
their Mellin transform and inverting the set of equations would lead to very complicated sets of 
equations when the weights become large. Hence the interesting functions are more or less the ones 
that have been built up from the original weight one functions by composing higher and higher 
integrals like F^ 0++0 (x) etc. without writing the result in terms of individual polylogarithms. 



7 Some applications 

The values at infinity of the previous section have some rather relevant applications for certain 
classes of integrals. This can best be illustrated with some examples. The following integral would 
under normal circumstances be rather difficult, but with all the above tools it becomes rather 
trivial: 

/ V n(x) ln2{1 ~ x) In(1+X) = -f)J- f'dx^Hx) ln(l-x) ln(l+x) (63) 
Jo x J^i+l Jo 

The integral is just one of the Mellin transformations, and hence the program will handle it. The 
sum is of the same type as all other sums in the Mellin transformation and hence will be done also 
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by the program. In the end the answer is expressed in terms of 5-functions at infinity which are 
maximally of weight 5 and hence they can be substituted from the tables. The final result is: 

1 , ln(x) ln 2 (l-x) hx(l+x) 3, , 2 o , . 7, , , 7 A 

dx ^J ^ 1 V )_ = _ CaCs _ C2ln 3 (2)+ C3 ln 2 (2)- -Cs 

X o o 4 I 

+41n(2) Li 4 (l/2) + ^ ln 5 (2) + 4 Li 5 (1/2) (64) 
15 

Similarly one obtains 

■ MX) Ib'(I-x) h'(l + x) = _i C2ln 4 (2) _12£ ( 3 + Z C3ln 3 (2) 

x 2 s w 140 s2 6 s v ; 



o 



37 Q1 

-^gCl - yC 5 ln(2) +81n(2) Li 5 (l/2) + 41n 2 (2) Li 4 (l/2) 

+1 in 6 (2) + 8 Li 6 (l/2) + 2 S^-i(oo) (65) 
y 



and the even more difficult integral 



/•i ln(l-x)Li 2 (^)Li 3 (^) 7 _ 5673 17 

I dx ~ ^ — = -^C2C3ln 2 (2)-^C2C5-5C2ln(2)Li 4 (l/2)- I ^C 2 ln 5 (2) 

1517 5 1 

-5C 2 Li 5 (l/2) + ^ClCa + gClln 3 (2) - -£]n(2) 

-^ C 3ln 4 (2) - ^ 3 Li 4 (l/2) - ^Cfln(2) - |c 5 ln 2 (2) 

74415 43 
+ T792" C7 " 181n(2) Li6(V2) " 14 ln(2) ^.-iM 

-61n 2 (2) Li 5 (l/2) -ln 3 (2) Li 4 (l/2) - ^-ln 7 (2) 

-24 Li 7 (l/2) - ySU, M (oo) + y 55^1(00) (66) 

As one can see, this technique allows the evaluation of whole classes of integrals that go considerably 
beyond the integrals in ref |fLlfl . 

Another application of the techniques of the previous sections concerns the evaluation of cer- 
tain classes of Feynman diagrams. When one tries to evaluate moments of structure functions in 
perturbative QCD one has Feynman diagrams which contain the momenta P and Q. Assuming 
that the partons are massless one has that P 2 = and because all dimensions are pulled out of 
the integral in the form of powers of Q 2 , there is only a single dimensionless kinematic variable left 
which is x = 2P-Q/Q 2 . The power series expansion in terms of P before integration corresponds to 
the expansion in terms of Mellin moments of the complete function after integration. The complete 
functions have been calculated for the two loop level |l5| but for the three loop level the calculation 



could only be done for a small number of fixed moments 2,4,6,8 and in one case also 10 [16]. To 
evaluate all these moments requires that the expansion in P should be in terms of a symbolic power 
N. This will introduce sums and these sums will be expressed in terms of harmonic series. After all 
integrals have been done all attention has to be focussed on the summations and it is actually for 
this purpose that the program SUMMER has been developed. By now a general two loop program 



has been constructed [17 and studies are on their way to create a three loop program. It should be 
noted that in the two loop program no series at infinity can occur. This puts a restriction on the 
functions that can occur in x-space. They have to appear in such linear combinations that all the 
constants (with the exception of £3 which comes from expansions of the T-function) should cancel 
in the Mellin transform. 
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8 Conclusions 



The algorithms presented in this paper provide a base for working with the sums that can occur in 
many types of calculations, one of which is the evaluation of Feynman diagrams in deep inelastic 
scattering. Additionally they allow the analytic evaluation of whole classes of integrals. The 
problem of the Mellin transforms of whole categories of functions has been solved, and a numerical 
solution for inverse Mellin transforms has been given. Most of the algorithms and tables have been 
programmed in the language of FORM version 3 and are available from the homepage of the author 
( http: / / norma, nikhef .nl/^t 68 / summer ) . 

The author wishes to thank D.A. Broadhurst, T. van Ritbergen and F.J. Yndurain for discus- 
sions and support during the various phases of this project. He is also indebted to S.A. Larin for 
the suggestion to have a look at these sums. 



A Conversion to the Basis 



To convert products of ^-functions with an identical last argument to the basis of single higher 
^-functions one can use a recursion. If one starts with the functions and and accumulates 
the results into the function the recursion reads: 



^Llj 1 ---j r ( n )Smlp 1 ---p 3 (n)Sq- l }..q t (n) — > S^^ 1 ...j r (n)Spi..^ s (n)SqJ.. qtm2 (n) 

- •> ^m 2 pi ■■■p s { n ) Sqi -q t rm{n) 

-^.^(nJ^WS^^^^ (67) 
The recursion starts with S^(n) = 1 and the recursion terminates when either S'W(n) or 

5(2) ( 

n 

has no more indices and hence can be replaced by 1 after which 

s£l jr (n)Sail qt (n) -> S qi ... qtjl ... jr (n) (68) 

with a = 1,2. Because this is a direct construction of the result, it is rather fast. It can be 
implemented in the language of FORM (version 3 or higher) very efficiently: 

repeat ; 

id,once,S(R(?a) ,n?)*S(R(?b) ,n?) = SS (R(?a) ,R,R(?b) ,n) ; 
repeat id SS(R(ml?,?a) ,R(?b) ,R(m2?,?c) ,n?) = 

+SS(R(ml,?a) ,R(?b,m2) ,R(?c) ,n) 

+SS(R(?a) ,R(?b,ml) ,R(m2,?c) ,n) 

-SS(R(?a) ,R(?b,ml*sig_(m2)+m2*sig_(ml)) ,R(?c) ,n) ; 
id,SS(R(?a) ,R(?b) ,R(?c) ,n?) = S(R(?b,?a,?c) ,n) ; 
endrepeat ; 

Note that the function SS carries the indices of and in this order. The function sig_ 

returns the sign of its argument. Hence the expression that uses this function is one way of writing 
the pseudo addition &. 

The above code has been made into a FORM procedure. A rather nontrivial test program could 

be: 
#- 

#include nndecl.h 
.global 
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L F = S(R(l,l,l,l,l),n)*S(R(-l,-l,-l,-l,-l),n); 
#call basis (S) ; 
. end 

It gives the result 

Time = 0.64 sec Generated terms = 1683 

F Terms in output = 1683 

Bytes used = 85104 

The run was made on a Pentium Pro 200 chip running the NeXTStep operating system. As one 
can see, these expressions can become rather complicated. On the other hand, weight 10 functions 
are of course not trivial. It should be noted that it is relatively easy to test routines like the one 
above. One can try them out for any functions and any values of the argument and evaluate the 
corresponding harmonic series into a rational number and see that they are identical. 



B Conjugations 

For the conjugations one should consider only S-functions with positive indices. The conjugation 
is defined with the sum 

(/("))° = -£(-!)*(?)/(<) ( 69 ) 
i=i 

That this is a conjugation can be shown easily by applying it twice. This gives the original function. 
For the function / one can use S'-functions or the combination of an S-function and a negative 
power of the argument of the S-function as in Sj 1 ...j r (n)/n k . For these functions one has: 
Theorem: The conjugate function of an element of the natural basis with only positive indices is 
a single S-function of a lower weight with only positive indices, combined with enough negative 
powers of its argument to give the complete term the same extended weight as the original function. 
Proof: First look at the weights one and two: 

1/n (70) 
Si(n)/n (71) 
1/n 2 (72) 

They clearly fulfill the theorem. Then write 

(S mjl ... jr (n)f = \{S h ... jr {n)/n^) c (73) 

This identity can be obtained by writing the outermost sum and then exchanging it with the sum 
of the conjugation. Assume now that the theorem holds for all functions with a lower weight. 
There are two cases: m = 1 and m > 1. When m = 1 the problem has been reduced to the same 
problem of finding the conjugate but now for a function with a lower extended weight. Hence, 
if the theorem holds for all simpler functions it holds also at the current weight. For m > 1 the 
conjugate of Sj 1 ...j r (n)/n m ~ 1 must be a single harmonic function of weight m — 1. This can be seen 
when one realizes that for each extended weight there are as many functions with their 'proper' 
weight equal to this extended weight as with their 'proper' weight less than the extended weight. 
Hence the function must have a conjugate that is a single S-function. Together with the fact that 
two conjugations give the original function, and the fact that all S'-functions of a given weight are 
linearly independent this completes the proof of the theorem. 



(Si(n)) c = 
(S2(n)f = 
(Si,i(n))° = 
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Next is the derivation of an algorithm to find the conjugate of Sj 1 ...j r (n)/n m ~ 1 . One way would 
be to successively build up the algorithm by first deriving all conjugates up to a given weight. After 
that one can obtain the needed conjugates by reading the formulae backwards. This is not very 
elegant. For a more direct way one can define the concept of the associate function. 

s i-j» = E(Su..*(0) c ( 74 ) 

1=1 

Note that because X = Sj 1 ...j r is an element of the basis, {X{i)) c contains powers of 1/i and the 
sum gives again a single harmonic function of the same weight as X. It is rather easy to prove 
that (X A ) A = X. The task of finding the conjugate can now be reduced to the task of finding the 
associate function. If this associate function can be written as S m j 1 ...j r (n) the conjugate will be 
Sj 1 ...j r (n)/n m . Similarly a function in combination with negative powers of n can be rewritten as 
a sum (Sj 1 ...j r (n)/n m — > S m j 1 ...j r (n)), and then the associate function of this function will be the 
needed conjugate function. 

The associate function can be found by construction. Assume that (Sj 1 ...j r (n)) A = S mpi ... Ps (n). 
Then 

n 

C 



(Si n - 3r (n)) A = £(S yi ... Jr (*)) 
i=i 

n -, 



i=i 1 



1 S Pl ... Ps (i) 



i 1 1 



— ^{m\-l)pi-p a { n ) (75) 



and similarly for k > 1: 

n 

(Skn-.. jr (n)) A = J2(S kjl ... jr (i)f 
i=i 

A 1 S h ... jr {n){i) c 

2_, i y ik-i > 
i=i 

n ^ 

i=l 1 

= S lqi ... qt {n) (76) 

with (S^_iy 1 ...j r (n)) A = S qi ... qt (n). Considering that (Si(n)) A = S\(n) associate functions to any 
weight can now be constructed. This algorithm is also easy to implement in a program like FORM. 

C Sums involving n — i 

In this appendix sums of the type 
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will be considered. It is impossible to combine the sums to a single basis element. Hence a different 
method is called for. Assume first that k > and m > (below). Writing out the outermost sum 
of the S- function with the argument n — i leads to 

E S mpi :.p s (n — i)S qi ... qt (i) y^y^ S Pl ... Ps (j) S qi ... qt (i) 

i=l i=l 3=1 •* 

y^ y^ S Pl ... Pa {j —i) S qv .. qt (i) 

y^ y~^ ^pi--ps(j~^) S qi --- qt {i) 
~ (j-i) m i k 

3=1 i=l w ' 
n i— 1 y 
= EE S Pl ... Ps (i-j)S qi ... qt (j) 7T—r— k (78) 
i=l j=l ^ -'J J 

Partial fractioning of the denominators gives sums in which the denominator is a power k' < k of 
j and sums in which the denominator is a power ml < m of i — j. These last sums can be done 
immediately by reverting the direction of summation. Hence: 

E Smp 1 -p a (n—i)S qv .. qt (i) _ ( k+m — 1 — a, „ f„-_„-w f-n _ 

jfc _ m _l ^ Zv Z^ ^Pi-PA' J) D qi-gtU) im ^-a - a 

i=l a=l i=l j=l ■> 

m , -. n i—1 1 

+ E( Ti )EEV- ft (j)V*(i-te(ra) 

a=l i=l 3=1 J 

Now the innermost sum is of a simpler type. Hence eventually one can do this sum, and after that 
all remaining sums are rather simple. Therefore this defines a useful recursion. When a negative 
value of mora factor (— l) 1 is involved things are only marginally more complicated. This algorithm 
has been programmed in FORM and carries the name sumnmii. An example of its application is 

#- 

#include nndecl.h 
. global 

L F = sum(j,l,n-l)*S(R(l,2,l),n-j)*S(R(-2,-l,-2),j)/j~2; 
#call sumnmii () 
. end 

with the output 

Time = 0.28 sec Generated terms = 478 

F Terms in output = 208 

Bytes used = 9148 

which are all terms with a single function of weight 11. 
The algorithm for doing the sums of the type 

= -^(-l)\ n .) Spi "^ n ~^ Sqi "^ (80) 
i=i 
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is more complicated. First one has to assume that all pj and qj are positive. Assuming also that 
k > 0, one can derive 

Gmp^~p s n ) = ^mpr'-ps 71— 1) + —Gmp^-ps (& — 1, n) 

I 1 v~y ntfi/ n \ S Pl ... Ps {n — i)S qi ... qt {i) 

n ^ ' [ i } ' (n-i) m -i k "' " y ' 

i=i K 1 

Because the weight of the G-function in the second term is one less, and because one can partial 
fraction the last term in the end all terms have a sum over a combination with a lower weight. This 
means that one can use this equation for a recursion, provided one knows how to deal with the case 
k = which is not handled by the above equation. For k = one obtains after some algebra 

rrm^qi-qt fn „\ _ i I W 1 ^ I U \ Sm lPl ---p s ( n — i)S qi --- qt (i) 

i=l 

i ^ W 1 \ffl / ' n \ 'S'pi---p s ( n ~^)'S'rra2gi---9t(^) (on\ 

i=i v ' 

Hence also here the weight has been decreased and one can use it for a recursion. The final 
expression for G can be obtained by an extra sum, because G(k, 0) = for all indices and one 
obtains an expression for G(k,n) — G(k,n — 1). One should also realize that in some cases it is 
necessary to change the direction of the sum {i — > n — i) which will introduce terms of the type 
(— l) n and hence this last sum can give S'-functions with a negative index. 

The routine that implements these algorithms (sumnmic) is a bit lengthy. A test run gives 

#- 

#include nndecl.h 
.global 

L F = sum(j ,l,n)*sign(j)*bino(n, j ) *S (R(l ,2, 1) ,n-j)*S(R(2,l,2) , j)/j~2; 
#call sumnmic () 
. end 

with the result 

Time = 0.36 sec Generated terms = 238 

F Terms in output = 131 

Bytes used = 6478 

and a simpler example gives 
#- 

#include nndecl.h 
.global 

L F = sum(j ,l,n)*sign(j)*bino(n, j)*S(R(2) ,n-j)/j ; 
#call sumnmic () 
Print ; 
. end 

F = 

- S(R(-3),n) - 2*S(R(-2,1) ,n) - S(R(l,2),n) - S(R(2,l),n) 

- S(R(3),n); 
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D Some sums to infinity 

There are special classes of sums for which the upper bound is infinity. A number of them can be 
evaluated to any level of complexity. Consider for instance the following sum (with p\ > 0; negative 
values just give extra powers of —1): 

F(m) = £ g *-*(>' + m ) (83) 
j=l 3 

Such sums can be evaluated by setting up a sum over m: 

° Sp 2 ...p r (J+m) 



F(m) = F(m-1) + J2 



pi 



= F(m-l) +Y,( Pl+ k_i l )^^( S i,P2,---,Pr(oo)-S ijP2 ,.. :Pr (m)) 
i=i 

* px+fc-l-t (-1)^ ^ Sp2 ... Pr (j+m) 

+ pi _i ) mPl ^ 2. * ( 84 ) 

i=i jf=i 

The sum in the last term is of the same type as the original sum, but it is of a simpler nature. 
The sum over i can just be worked out, because k and pi are just numbers. Hence this defines a 
recursion which can be worked out, if not by hand, then by computer. In the end one obtains an 
expression for F(m) — F(m — 1) which can be summed: 



F(m) = F(0) - F(i-l)) 

i=i 

m 

= S kpi .. Vr (oo)+'E(F(i)-F(i-l)) (85) 



i=i 

Similarly one can consider sums of the type 



oo 



The technique to construct a recursive solution for these sums is similar. One can study the function 

F(m) - F(m-l) = V 5 Pi-Pr(j) _ y ■V-^O' + l) 

p[ (j + m) k p Q (j + m) k 



S Pl -Pr(l) ^ ( Pl + k-l-i (-l)Pl-^ 

,* 2^ fc-1 J( m _l)PiH^ Di iP2>-iP'- l > 00 ' 1 

j=l ^ ' 



^ Pi- 1 (m-l)P!+^ ^ (j+m)* 1 j 

and again the last term is of a simpler nature. Hence there is a useful recursion and these sums 
can be solved. 

In both cases there will be some S'-functions in the answer that have the argument infinity. 
These should not present any special problems as they have been discussed before. 
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E Miscellaneous Sums 



In this section some sums are given that can be worked out to any level of complexity, but they 
are not representing whole classes. Neither is there any proof for the algorithms. The algorithms 
presented have just been checked up to some rather large values of the parameters. 

The sums that are treated here involve two binomial coefficients. There are quite a few of these 
sums in appendix ||], but here are the ones that can be done to any order. The first relation that 
is needed is: 

jLHn l+j X m+i+j ) = (-1) ( i )( i ) 

Taking m = n — i leads to: 

n . n ■ 

for < m < n. Here / indicates the conjugation of appendix |B[ This is a rather useful identity 
as it divides the necessary amount of work by two. Alternatively it may even make terms cancel 
and hence make further evaluation unnecessary. 

A new function is needed to keep the notation short: 

U k (n,m) = S k (n+m) - (-\) k S k (n-m) - S k (m-l) (90) 

for k,n > and m > 0. Uo(n, m) is defined to be one. 

One of the ways the harmonic series can be introduced in many calculations is by expansion of 
the T-function. At the negative side its expansion is: 

T(-n + e) = ^-r(l + e)(l + 5 1 (n)e + 5 1 , 1 (n)e 2 + 5 1 , 1 , 1 (n)e 3 + 5 1 , 1 , 1 . 1 (n)e 4 + ---) (91) 
en! 

Actually these special harmonic series can be written as a sum of terms that contain only products 
of harmonic series with a single sum as in: 

25i,i(n) = (Si(n)) 2 + S 2 (n) (92) 
65i,i,i(n) = (S 1 (n)) 3 + 35i(n)S 2 (n) + 2 < S 3 (n) (93) 
24S r ll i,i ) i(n) = (S 1 (n)) 4 + 6(5 1 (n)) 2 5 2 (n)+85 1 (n)S 3 (n) + 3(5 2 (n)) 2 + 6S4(n) (94) 

Notice that the factors are related to the cycle structure of the permutation group. One can define 
the higher U functions by analogy: 

2J7 lia (n,m) = {U x {n,m)) 2 + U 2 (n, m) (95) 

6^1,1,1 (n,m) = (Ui{n, m)f + 3Z7i (n,m)U 2 (n,m) + 2U 3 (n,m) (96) 
24L r i i i i i i i (n,m) = (Ui(n, m)) 4 + 6(Ui(n, m)) 2 U 2 (n, m) + SU^n, m)U 3 (n, m) 

+3(f7 2 (n,m)) 2 + 6C/ 4 (n,m) (97) 

With these definitions one can write (0 < m < n): 

^ ^ fl — |— j \ fl I ( fYl — X ^ ^ 

f^}~ 1 ^^ 3^m+j\m+jY = (n + m)! (98) 
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t(-iyC)( n+J .)— 1 —, = n \ (m ~!, )! U^m) (99) 

±(-iyC)( n+J .)— ^3 = "1 (m ~!, )! £/n(n,m) (100) 

^ J rn+j (m+j) 3 (n + m)! ' v ' ; v ; 

EW( n )( n+J )^4 = "! (m "!, )! ^in(n,m) (101) 

^ J m+j (m+j) 4 (n + m)! ' ' v ; y J 

etc. In the case that m is zero there are different expressions: 

n . -, 

EwCX;)- = - 2S i(») ( 102 ) 

E(- 1 ) J (p( n 7 J )-2 = -4Si,i(n) + 2S 2 (n) (103) 

n . -. 

EHXCX ; )-3 = -85i,i,i(n)+4Si >2 (n)+4S 2> i(n)-253(n) (104) 
i=i J J J 



n . -■ 

SHAN'T )-4 = -1651,1,1,1 (n) + 8(5i,i, 2 (n) + S 1>2>1 (n) + S 2 ,i,i(n)) 

-4(5i, 3 (n) + S 2 , 2 (n) + 5 3>1 (n)) + 2S 4 (n) (105) 



j=1 J J J 



and the pattern should be clear: For l/j k there will be all functions with weight k. The ones with 
m nested sums have a coefficient _(_i) fc - Jm 2 m . A recipe of a similar type is found for the following 
sums: 

n . 

EWCX' j)Si(j) = 2(-l)»Si(n) (106) 
j=i J 

n 



,n.,n+j 

3=1 

n 



E(-1) J G)( j ) S 3(J) = (-l)"(25^(n)-45_ 2 ,i(n)) (108) 

n . 

EHX'CX; )^4(J) = (-l)"(-2^(n)+4(^,i(n)+&. 2 , 2 (n))-8S_ 2 ,i,i(n)) (109) 
j=i J 

and the pattern here is that one should make all higher series that start with a negative index that 
has a value of at most —2, after which there are only positive indices. All functions are of weight 
k (for Sk inside the sum), and for m nested sums the coefficient is (_ i )"+fc- jm 2 Tn . The exception is 
k = 1 but that is because S\ is its own associated function and its conjugate is purely of the type 

i/jK 



F Summation tables 

During the work that inspired this paper quite a few other sums were evaluated that are not 
represented by the above algorithms. Many of these sums can only be done for a fixed weight and 
most of them were not readily available in the literature. Hence they are presented here in a number 
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of tables, even though eventually many of these sums were not needed in the final version of the 
program. For completeness also a large number of sums are presented that are already available in 
the literature. A number of these sums can be derived formally. Some were derived by 'guessing' 
and then trying the resulting formula for a large number of values. 

In all sums it is assumed that all parameters i,j,k,l,m,n are integers and have values > 0. In 
some cases the formulae can be extended to noninteger values. 

It should be noted that all sums that can be handled by the procedures of the previous appen- 
dices are not in the tables. They would make the tables unnecessarily lengthy. 

Some formulae that are used very often are presented first: 

£(m+i)l = (n+m+1)! 
f-f i\ n\(m+l) 

E{m+i)\ _ (n+m+1)! ml . . 

_ (k+i)\ ~ (n+fc)!(m+l-fc) ~ (fe-l)!(m+l-A;) ^ ' 



= S{n) (112) 

3=0 J 

EH^") = H) m ( n ~ ) ( 113 ) 

3=0 J 

±C.){m+j)\ {k+n-j)\ = rn\M^^±^ (114) 



3=0 



v ■ n (m+j)\ (w+fc-1)! m! 

^ K ' { j'(m+k+j)\ (fe-1)! {m+n+k)\ 1 ' 



3=0 



The last three formulae can be extended to noninteger values of m and k. They can be used 
occasionally before T-functions are expanded to yield harmonic series. 
At times some auxiliary functions were needed. They are defined by 

n -I 

A k (n,m) = ^(-l)^( n + m ) (117) 
j=i J 

b a+b—i 

A(l,l,a,6) = EE" ( 118 ) 

i j 

i=a j=a u 

Rm,k(n) = (119) 

Sometimes A is not always the easiest function to manipulate. Therefore the function A' is some- 
times handy: 

A'(l,l,a,6) = ^(6-a)(A(l,l,a,6)-(5 1 (a-l)-,S 1 (6)) 2 ) 

-i#(a-6-l)A(l,l,6+l,a-l) (120) 

Both functions involve summations over triangles in the two dimensional plane. These triangles do 
not touch the origin. 
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F.l Sums without (— l) j 

First is a number of expressions that are at the lowest level of complexity. 

Ai(n,l) = -Sifn+lJ + tff (121) 

Mn,2) = - Slin+2)+ ( Jjr +2) +Hr ,122) 

Ai(l,m) = -(m+1) (123) 

Ai(2,m) = -(m+2)(m-3) (124) 

Ai(n,m + 1) = ^(n,m) + (-l) n ( n+m ) 3-— (125) 

Ai(n + l,m) = ^i(n,m + l)-(-l) n ( n+m+1 )^ (126) 

m n+1 



A k (n,m+1) = J2(n+m+iy- k Mn, m) + ((-l)"( n + m ) - 1) (n+ ^ +1)fc (127) 



A fc (n,m+1) = A fe (n,m) + ^(-iy( n+m+1 )S jfc (j)-(-ir( n+m ) 1 S fe (n) (128) 

j= i 3 

Next are some miscellaneous sums: 

E-^i(6+j) = (Si(a)-5i(6))S , i(a+6) + 25i i i(6)-^5 2 (6) 

3=1 3 

+A'(l,l,a,6) (129) 

E( n "! +J )^T7 = ( n+ 7i m )^--(-l) n A 1 (n,m)-(-irS 1 (n+m) (130) 
J n-\-j n n 

±^±A = 2S 1>1 (n)-±S 2 (n) (131) 
^^J) = lS - 1 ( 2n )5 1 (n-l)-25i ) i(n-l) + ^5 2 (n) (132) 
Y, S - l(n , +j) = S-i,-i(n) + Si,-i(n)-ls 2 (n) (133) 

j=l ^ * 

E 5l(n 2 +j) = 5i, 2 (n)-25 2) i(n)-^ 3 (n) + 212 1)2 (n) (134) 
i=i 3 

jr S - l{n 2 +3) = 2S 2)1 (n)-fiU 2i _ 1 (n) + &. 1 ^(n)-52 1 -i(n)-^5 3 (n)-2i2 1)2 (nXl35) 

E g " 2(re+J) = ^ 2 ,_ 1 (n) + 5 li _ 2 (n)+5 2 ,_ 1 (n)-5 2il (n)-^ 3 (n)+ J R 1 , 2 (n) (136) 
m i m -i 

]T — — -Si(j) = S l (m)Si(n+m+2)-Y / Si(n+j + l) (137) 
j=l n+] + i j=1 J 
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m -i m -i 

E— + = - £ — -^(n+l+j) + < S 1 (m+l)5 1 (n+m+2) - (138) 

The last two equations are not solving anything, but they are useful in the derivation of some of 
the next sums. First an equation that is like a partial integration. 

£( m + J )/(n+W) = £( m ? +J )E/« ( 139 ) 
j=o J j=0 J i=l 



It is central in the derivation of the next equations 

i+js 1 ,n+m 
j n+l-j n + 



tC+ j )^- = (^^(n+m+l)-^)) (140) 
j=0 

+2 S 2 (n+m+l) - 2S 1 , 1 (n+m+l) + ]T Si{n+m+l-i) ) (m) 

i=i 1 

U )Sk{n+l-j) = )5 fe (n+l-j) 
j=o J m+i j=0 J 

E( ro 1" J ')Si(n+l-j) = ( n +^+ 2 )(5 1 (n+m+2)-5 1 (m+l)) (143) 
j=o 3 n+i 

^(^^^(n+l-i) = I( n +™+ 2 )(2Si,i(ro+l)- 25 1 , 1 (n+l) + S 2 (n+1) 



(n + m+2)(m+l) 



A'(l,l,m,n+1) - A'(l,l,m+l,n+l) 



+ (25i (n+m+2) - _l--)( < S 1 (n+l) - 5i(m+l)) ) (144) 



£( m + J )<? 2 (n + W) = ( n + ^ + 2 )(5 2 (n + m + 2)-5 1 , 1 (n+m + 2) 



+S 1 (n+m+l)S 1 (n+l) + gl(m) n + ^S 2 (n+1) 

n+m+2 2 

-5i,i(n+l) - A'(l,l,m,n+1) ) (145) 



Similarly one can derive: 



YS m ^ J )Si(m+j) = ( n+ ^ +1 )( 5l (n+m+l)-^) (146) 

V( m+J )^(m +J ) = ( n+TO+1 )((S 1 („+ fn) - J_) 
J n yu iV ' m+l n n+m+l m+l ; 

+5i i i(n+m))-(-l) ro ^- I (Si(n+m) + Ai(m,n)) (147) 
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Y J { m+3 )S 2 (m+j) = ( n+m+1 )S 2 (n+m)- { -^(S 1 (n+m) + A 1 (m,n)) (148) 
~^ j n m+1 

These formulae give also some 'partial integration': 

^(^l Sl (m+j)f(n-j) = g m (m+i-l)! ^ (m+j _ 1) g /( . ) 
j=o j=o i=l 



j! m 



E" li nh S 1 (m+j)f(n+l-j) = -( m -l)J2 { - ry±^Si(m-l+j)/(n+2-i) 



3=0 J ' j=0 J 

n4-l 



+ £" ^Sx(m-l+i)/(n+2-i) (150) 



7'! 

i=o J 



This is used for the next equations 



g ( m+i^(rn+i) = ( m+n _ 

+5 2 (m) - 5i(m+n)5i(m)) (151) 

2( m + J )5i(m+i)^(n-i) = ( n + ™ +1 )(-(S 2 (n + m+l)-S 2 (m+l)) 
i=o 7 

+ (5i(n + m+l)-5 1 (m+l))(5i(n+m + l) ?—)) (152) 

m+1 

F.2 Sums with 

The next sums all contain a factor (— 1) J and hence they give completely different results than the 
corresponding set of sums without the (— 1) J . The most important ones have been treated in the 
appendices || and |C|. 

V( n )(-iy-^ = ( n+m )-i± (153) 

~q 3 m+j n ' m y ' 

f2( n )(-iyS 1 (m+j) = -( n + m )- x - (154) 
i=o J n n 

E(pH) J (^-)f = C^r^iSAm+^-S^m- 1)) (155) 

= ( n + m ) -ii( 5l(m+n) -5 l(n )) (156) 
J m+j n m 

J2( n .)(-iyS 2 (m+j) = -{ n+m )-^hs l (n+m)-S l (m)) (157) 

3=0 J 

n -, 

^(^(-ly^iCm+i) = -( n + m )-i-(5 1 (n+m)-5 1 (n-l)) (158) 
i=0 J n n 
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EHH")^ = -A 1 (n ! m)((n + m)( n+ ^- 1 ))- 1 (159) 

There is also a number of sums with more than one binomial coefficient. First a few general 
ones. 

n ■ 

i=o J 2 

n . 

EH) i (")( = m<n (161) 

i=o J J 

n , ■ . -I 

E(-iy'(")( m | J )j = (-<>™>™-i 

= m<n-l (162) 

n . ■ 

EOC+jhx = 0<m<n ( 163 ) 

i=o 7 ^ 

ECx^+^x-iy = H) n ^>o (164) 

^j!(n+l-j)! (n+l-j)! 

V— — - (n ~ 1+ { )! (-iy = -6(n-l)--5(n-2) (166) 

^il(n-l-i)! j!(j + 2) y ! 2 y J 6 K J y 1 

And now the sums with two binomial coefficients and an argument j in the extra piece (see also 
appendix |E|): 

n . 

£0( n 7)(-l)^i,i(i) = (-l) B (4S ltl (n)-2S 2 (n)) (167) 
E(p( n t J X-l) J 7 S iO') = 2S -2(n) (168) 



,n, ,n+j 

\ 

3=1 

n 

,n.,n+j 



Efyij Si,2(j) = (-l) n (2&_ 3 (n)-45 l5 _ 2 (n)) (169) 



E( '•)( 7 ^W^iO') = ("l) n (2^N) (170) 



E( ? J)( n t i )(- 1 )^i,i,i(j) = (-l) n (25 3 (n)- 4Si, 2 (n)-4S 2 ,i(n) H-S^n)) (171) 
i=i J 



n . -. 

E(!)( n 7 J )(- 1 ) J - 5 2(i) = -25 8 (n) (172) 



n . ■ -I 

EO( ^X- 1 ) V 1 ' 1 ^ = 4^ 2il (n)-25_3(n) (173) 
i=i J J J 

n . 

E(^)( n 7 J )(-l)^ S i0') = 45 lra (n)-2a fl (n) (174) 
i=i J 7 
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Similarly there are formulae with j + 1 (see also appendix ^). 
n ■ i 

^ ( i )( i+i )( 1} ^TT- " - ( - 1} ^TIF 

E(")(^)(-l)^i,l(i+l) = ^(Si(n+1) + Si(n-1)) 

n . -. 

D")(" + 1)H)'^ - -*(-«) -*(»-!)) 



J'=0 



1 1 

+ 



, =0 'J i + 1 i + l ra 2 (n + l) 2 (n + 1) 3 

E(p(^)H)^3(j + l) = ^^(l-SxCn+l)-^^-!)) 

And the formulae with j + 2 (see also appendix 

x(5i(n+2) + 5i(n-2))(£Li(n+2) - S-i(n-2)) ) 



j=0 



" ,n w n+j w 1 s,S 2 (j+2) 1 



5 1 1111 5 1 

+ - 77 7 + 



The formulae with j + n: 



6n-l 2n 2n+l 6n+2' 



V( n )( n+J )(-ir^- = o 

p 3 J n+j 

n ■ 

^Qn J )(-l) j Si(n+j) = 2(-l)»S 1 (n) 

A .n u n+j v v - 1 , (re-l)!(n-l)! 

J )H) (^7F = " H) (2n)! 
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f2( n ,)( n ^ j )(-l) j (S 1 , 1 (n+j)-S 2 (n+j)) = (-l) B (4Si,i(n)-35 2 (n)) (190) 

3=0 J J 
n . 

£(")( i 3 )HY-Mn+j) = 35 2 (n) + 25_ 2 (n)-45 1)1 (n) (191) 

Djx^JHy^acn+i) - (-D"^i) (192> 



j"=o 

n 



? ( j )< i+i )H),Sl( " +j) = H) "STT (193> 

g ( » )( - j)( _ 1)JS! , !+J) _ 



J=0 



1112 1, , N 



?'i )( "i i)H)i (FiF Sl( " +i) " ^TI) 5 ^"- 1 '-'- 1 '"^? (197> 



n 



J n ' j + 2 v J> v 7 v 6n-l 2n 2n+l 6n+2 7 v ; 



i=o 



B^Ohi^^') - (^< s '<"- 2 >+ 2 H)> 2 + » + l)^§) < 199 > 



i=o 



E(p(3+ J J )H) J ^i(^) = ^j(5i(n-3)-(SLi(n+3)-flLi(n-3))) (200) 



3=0 



n . -I i 

E(")(i+j ^1)2 = ^ T (25 1)1 (n-l)-25 2 (n-l)-5_ 2 (n+l) 



j"=o 



-6L 2 (n-l) + gi(n+l)5i(n-l)- ( ^ - ) ) (201) 



^ )( 2+j )(_1 ^ Sl(n+i) = (^yy( 2 ^i(— 2 )- 2 ^(— 2)-^i(n-2) 



-SL 2 (n+2) - SL 2 (n-2) + 5i(n+2)5i(n-2) 



+ (_l)n (n 2 +n+3) ^_|| } (2Q2) 



And some mixed formulae. 



V( n )( n+J )(-iy— ^- = (m<n) 
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g(p(T^(^ = _ Hr (mzl)!^^ (f| , ra) (204) 

E(p(^)(-l)^ 2 (m+i) = (_i r+ ^ (^-^ ! ( g_ l(n+m) _ g_ l(n _ m)) 

(n > m) (205) 

E(")(^)(-l)^i(n+i) = -H) ^^) 1 ! 1 (n " m) (2 ° 6) 
There are also sums with two factorials in the numerator. 

±ji(m-m-iy = ^±^ + (-ir (n+1)! (m ~ n)! (207) 

v ^ v y m+2 v ; m+2 



i=o 



m+2 v v ' m+2' 

(m+1)! 
(m+2) 



2 



(208) 



(n-j)\(-iy-^- = -(n+1)! (S 2 (n+l) + 2&_ 2 (n+l)) (209) 

£j! (n-i)!(-l)^ = n! (S 2 (n) + 2S_ 2 (n)) (210) 

E(-l) j i' ( ot m 7 )! = ("l) m (^+2)! (n-2)! (^(m+3,n-2) + S 1 (m+2)) 
j=0 ^ + m ^ 



+ ( ^3 )2 )! + (-l) m (-+l)! (n-1)! (211) 



^ ( _l ) ^ (n+m _ J - )!i , AxCj.n-j) = (^+^+l)' ( g i(n+m+1) _ g l(m+1)) 



3=0 



F.3 Partial sums 



n+m+2 

(n+l)!m! 

v y n+m+2 w v ; 



For derivations the use of partial sums (sums of only part of the range of the binomial coefficients) 
are very useful. Unfortunately these are hard to obtain and hence only a limited number of them 
can be presented here. 

m -. 

EH) J '0 = H) m ( n " ) ( 213 ) 

m -. 

EH) J '(,)^ = ^Kn-m) (214) 

m — 1 1 1 

EH) 3 (>iO') = HH* )(Si(m) + -)-- (215) 
m — 1 1 

^(-l) j OSi(n-j) = (-l) m C i )(S 1 (n-m-l) + -) (216) 
j=o 3 m n 
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